A new effective interaction for the trapped Fermi gas 



O 

o 

(N 

G 
3 



Y. Alhassid, 1 G.F. Bertsch, 2 and L. Fang 1 
1 Center for Theoretical Physics, Sloane Physics Laboratory, 
Yale University, New Haven, CT 06520 
2 Department of Physics and Institute of Nuclear Theory, Box 351560 
University of Washington, Seattle, WA 98915 
(Dated: June 27, 2007) 

We apply the configuration-interaction method to calculate the spectra of two-component Fermi 
systems in a harmonic trap, studying the convergence of the method at the unitary interaction limit. 
We find that for a fixed regularization of the two-body interaction the convergence is exponential 
or better in the truncation parameter of the many-body space. However, the conventional regu- 
larization is found to have poor convergence in the regularization parameter, with an error that 
scales as a low negative power of this parameter. We propose a new regularization of the two-body 
interaction that produces exponential convergence for systems of three and four particles. From the 
systematics, we estimate the ground-state energy of the four-particle system to be (5.05 ± 0.024)fttj. 
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PACS numbers: 31.15.-p, 03.75.Ss, 21.60.Cs, 71.15.Nc 

The study of cold trapped atomic condensates has be- 
come a rich field experimentally. By providing a strongly 
interacting system that is well defined, it also offers physi- 
cists an unprecedented opportunity to assess theoretical 
techniques that cross the boundaries of disciplines. In the 
so-called unitary limit, the only dimensional scale of the 
problem is fixed by the harmonic trap frequency. System- 
atic studies have begun on small systems using fixed-node 
Monte Carlo [l|, l2fl and density functional methods [3j. 
Remarkably, the exact wave functions and energies of the 
A = 3 system are known, calculated by solving a single 
transcendental equation 4] . Our work here is in the con- 
text of the configuration-interaction (CI) method, widely 
used in atomic, molecular, and nuclear spectroscopy. We 
study the convergence of the CI method with respect to 
a regularization parameter of the two-body interaction 
and find that a simple regularization scheme that renor- 
malizes the interaction produces slow convergence of the 
three- and four-particle spectra. We introduce a new ef- 
fective interaction that gives exponential convergence, at 
least in small systems. 

Hamiltonian. The cold trapped atom system is modeled 
by the Hamiltonian 
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where A is the number of atoms, to is the trap fre- 
quency, and Vq is the interaction strength. We have two- 
component fermionic systems in mind, which controls the 
symmetry of the allowed states. The interaction is rep- 
resented as a S function (contact interaction) but as we 
shall see below it requires a regularization. Here we focus 
on an attractive contact interaction in the unitary limit 
of infinite scattering length. 

The two-particle problem. The two particle system (A — 
2) is separable in center of mass and relative coordinates 
r = T2 — i"i . The center of mass Hamiltonian describes an 



harmonic oscillator with frequency u> and mass 2m, while 
the relative-coordinate Hamiltonian is H Ic i = — |^V 2 + 

mj 2. The two- 



^liu> 2 r 2 + Vb<5(r) with reduced mass p. 



particle energies are given by 



E = {2J\f + £ + 3/2)huj + e n i 



(2) 



where A/", C and n, I are the radial quantum number and 
angular momentum of the center of mass and relative 
motion, respectively. The energies e n i are the eigenvalues 
of -f/rei, and may be derived from the boundary condition 
at the origin imposed by the unitary interaction [j| . The 
contact interaction affects only the I = partial waves, 
and shifts each s-wave oscillator energy down by one unit 
of fajj 6| . Thus we have 



a = (2n + I + 3/2 - 8i fi )frw 



0,1,2, 



(3) 



The renormalized contact interaction. In the CI method, 
the contact interaction in Eq. ([1]) must be treated explic- 
itly. However, a (5-function interaction cannot be used in 
three dimensions without a regularization. We shall do 
this by truncating the space of relative-coordinate wave 
functions to a q subspace defined by the lowest q + 1 os- 
cillator I = wave functions (see also Ref. 0)- Within 
the truncated space the relative-coordinate Hamiltonian 
can be written as 



(i? re i)i^ = (2n+3/2)fio;5 n ,„/ + yW, {0<n,n'<q) (4) 
where 

vS = fiwx,iMoM»'(o), (5) 



r(9) 



and V™(0) = 7r- 3 / 4 v /(2n+ l)!!/(2"n!) is the (n,l = 0) 
oscillator wave function at r = for an oscillator of radius 
1. The parameter \q is a dimensionless normalization 
constant related to V by \q = (ft 2 /Ai) _3/2 (ftw) 1/2 V . 

We determine the normalization constant Xq by requir- 
ing the ground-state energy of the truncated Hamiltonian 



2 



to equal the exact value for the unitary contact interac- 
tion, £00 = fvuj/2. The separable form of ||SJ| permits 
an algebraic diagonalization of the Hamiltonian. Each 
eigenvalue e of (UJ) satisfies the dispersion formula 
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Requiring e = eoo = ftw/2 in ([6]), we obtain a closed 
expression for the normalization constant 



Xq 



T 3/2 



(2n- 1)!! 



(7) 



We note that the sum in ([7]) diverges as q 1 ! 2 for large q Q. 
Thus, the strength of the ^-function goes to zero as q — > 
00, showing the need for a renormalization procedure. A 
similar relation between the strength of the interaction 
and the cutoff can be derived for a plane-wave basis. In 
that case the relation is Vq — —Tr 2 h 2 /iJ,A where A is a 
momentum cutoff |9j. This value of Vo agrees with the 
asymptotic expression of Eq. ([7]) || once we equate the 
corresponding cutoff energies as h 2 A 2 /2[x = (2g+3/2)fiu. 

The excited states of the g-truncated Hamiltonian (j4|) 
have energies e^Q that differ from the exact unitary spec- 
trum ([3]). Using the dispersion relation ©, we find that 



the error in the energy Ss n 
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Sno goes to zero at 



large q, but only at a rather slow rate, Se^Q ~ g -1 / 2 . We 
present evidence below that this slow convergence is also 
present in the g-renormalized energies for the A = 3 and 
A = 4 systems. This makes it problematic to extrapolate 
the q series to estimate the true q — ► 00 energies. 
A new effective interaction. We have considerably more 
freedom to construct the g-space interaction than we have 
exploited so far. The only requirement on the g-space 
Hamiltonian is that it converge to the unitary limit for 
large q. For example, in effective field theory one may 
introduce derivatives of the contact interaction to fit cer- 
tain properties of the two-particle Hamiltonian. Here we 
propose the following prescription to improve the g-space 
interaction: simply require that the relative-coordinate 
Hamiltonian reproduce all q + 1 s-wave eigenvalues of 
Eq. Q . We can do this and still keep the separable form 
for the interaction, 



v:y = -tout* 



(8) 



A motivation for preserving the separable form is given 
in the discussion below. There are q + 1 independent 
variables /„ in the interaction |8]) and the same number 
of eigenvalue equations having the form of Eq. ([6]) with 
/„ replacing ^/|x^|V'n(0)- Using the conditions that all 
q + 1 lowest I — unitary eigenvalues ^ (n = 0, . . . , q) 
are reproduced, we find the following q + 1 equations for 



q 

11 = V 



r) + l 



1 (r = 0,...,g) 



(9) 



Eqs. © determine a unique solution for f 2 (n — 
0, • ■ ■ , q) [nj. We choose the sign of the real numbers 
/„ to coincide with the sign of VVi(O). Using the con- 
vention that the harmonic oscillator wave functions be 
positive at the origin, the unique solution for /„ is 
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'(2n+l)!! [2(g-n)-l]! 
(2n)!! [2(9 -n)]!! 



(10) 



The interaction defined by |8j) and ([10]) is different from 
the renormalized contact interaction for any q. However, 
its eigenfunction components (in the 3-D oscillator ba- 
sis) converge to the corresponding unitary eigenfunction 
components in the limit of large g with an error of ~ g _1 . 
In comparison, the eigenvector components of the renor- 
malized contact interaction converge to the same unitary 
eigenvector components but at a slower rate of ~ g~ 1//2 . 
CI method and truncation of many-particle space. In 
the CI approach, one uses a single-particle basis in the 
laboratory frame and constructs a many-particle basis 
of Slater determinants for A fermions. In our prob- 
lem, a natural choice for the single-particle basis are 
the eigenstates of the three-dimensional harmonic oscil- 
lator. These states are labeled by orbital quantum num- 
bers a = (n a ,l a ), the orbital magnetic quantum number 
m a , and an additional two- valued quantum number (e.g. 
spin) to distinguish the two species of fermions. 

A way to truncate the many-particle space must be 
specified, because there is no natural truncation associ- 
ated with the interaction except in the trivial cases q = 
or A = 2. There are a number of truncation schemes 
in the literature; here we will define a truncated single- 
particle orbital basis and construct the A-particle wave 
function allowing all possible anti-symmctrized product 
states. In particular, we shall use all single-particle states 
in the oscillator shells N = 0, . . . , A max with N = 2n a +l a 
to construct the many-particle states. There will be 
two limiting processes necessary to calculate the many- 
particle energies. The first is N max — > 00, which we will 
investigate for fixed q. Then, with converged g-regulated 
energies we estimate the q — > 00 limit. 

Two technical aspects of our calculations should be 
mentioned. The two-particle matrix elements of the in- 
teraction in the oscillator basis are conveniently calcu- 
lated using the Talmi-Moshinsky brackets to transform to 
relative and center of mass coordinates 11( . The many- 



particle Hamiltonian is constructed and diagonalizcd us- 
ing the nuclear shell model code oxbash [12j. Unlike 
the nuclear shell model, our orbitals are characterized 
by integer angular momentum values. The two fermion 
species are distinguished in the same way as neutrons and 
protons are distinguished in the nuclear application. 




FIG. 1: Convergence in -/V max for the A — 3 ground-state 
energy, (a) Ej$ versus iV max for q — 3. Open circles corre- 
spond to the renormalized contact interaction and solid circles 
to the interaction defined by and (|10|) . (b) AE^ versus 
A^max in a logarithmic scale. All energies are in units of Hu>. 

A=3 system. We now show the results for A = 3. The 
ground state of the A = 3 system is a negative-parity 
state with total an gula r momentum L = 1 and energy 
4.2727243 ... tuu [4j, |l_3j. In our CI convergence studies 
we computed the ground-state energies Epf for q — 
1,2,3,4 and N max = q,...,7. 

For a fixed q, we find that converge exponen- 

tially or better in A/ max for both interactions. This is 
demonstrated in Fig. [TJ Panel (a) shows Effl versus 
-^max for q — 3. Both the renormalized contact interac- 
tion (open circles) and the new interaction (solid circles) 
are monotonically decreasing, as they must when the 
space gets larger. The important point, seen in Fig.[]Tb), 
is that the energy differences AE$ = E$ _ 1 — 
decrease rapidly on a logarithmic scale. In fact, the de- 
crease is steeper than linear on that scale, suggesting 
that the convergence might be faster than exponential. 
The solid lines are quadratic fits to log(Ai?j^ ), used 
to extrapolate to a value of E^ = E& . We observe 
the decrease rate of AEpf to be monotonically in- 
creasing with Af max , so a conservative lower bound in 
E^ is obtained using a fixed-rate extrapolation above 
A/max = 7 with an average rate determined by the points 
N max = 5,6,7. An upper bound for E^ is given by E^ . 

Fig.rjja) shows the converged or extrapolated energies 
versus q. These energies are monotonically increas- 
ing function of q. For the new interaction (solid circles), 
we observe a fast convergence to the known exact value 
(dotted line) . Fig. Hfb) shows the absolute value of the 
deviation SE^ = E^ - E^ from the exact result in 
a logarithmic scale. The concavity of the curve for the 
renormalized contact interaction (open circles) indicates 
the convergence in q is slower than exponential. We find 
this convergence to be consistent with a low negative 
power law ~ q~ a with a in the range ~ 0.5 — 1.5 (for 
the excited A — 2 system it can be shown analytically 



FIG. 2: Convergence of the (/-regulated energies for the A — 3 
ground state, (a) E fq ' versus q for both interactions (symbols 
and units as in Fig. [1} . The dotted line is the exact ground- 
state energy, (b) The error [SE 1 -^] in a logarithmic scale. 

that a = 1/2). However, for the new interaction (solid 
circles) the convergence is at least exponential. 

This exponential convergence allows for an accurate 
estimate of E^°°\ We calculated successive energy differ- 
ences AE^ = E( q -V - and determined an average 
rate of decrease A of \AE^\ for q below a given q' . As- 
suming a fixed rate A for q > q' , the extrapolated energy 
is [XE^ - E^'-^/iX - 1). We can take this value to 
be an upper bound for E^°°\ since the rate of decrease 
of |Ai?( 9 )| seems to be a monotonically non-decreasing 
function of q. Using q' — 3 and an average decrease rate 
of 3.28 (determined from AE^ at q = 1,2,3), we find 
E (oo) _ ( 4 2?4 ± .004)fto;, an accuracy of 0.1%. 

We carried out a similar study for the LP = + first 
excited state at E^ = 4.6662 ... hw 0, [ljjj]. Results 
are shown in Fig. 3(a). As in the ground-state case, we 
observe a low negative power law convergence for the 
renormalized interaction and exponential convergence for 
our interaction. Using E^ and an average decrease rate 
of 1.83 obtained from q = 1,2,3, we estimate E^°°' — 
(4.646 ± Q.025)huj, an accuracy better than 0.6%. 
A =4 system. We also studied the L — ground state of 
the A = 4 system with two particles of each species. The 
results for E^ are shown in Fig. [3^b). Here the exact 
value is unknown. An upper bound (using the new 

interaction) is E^ = 5.075 fku. A lower bound can be 
obtained as for the A = 3 system. The inset of Fig. |3](t)) 
shows AE^ in a logarithmic scale versus q for the new 
interaction. Again, the convergence seems to be at least 
exponential. The straight line is a fit to log(A£(«)) using 
q = 1, 2, 3, and provides an average decrease rate of 2.14. 
Using the extrapolated E^ = (5.074 ± 0.001)^ and 
this average rate, we estimate E^ — (5.051 ± 0.024)ftu;. 
Our result agrees with fixed-node Monte Carlo estimates 
of (5.1 + ±0.1)7ttJ [l| and (5.069 ± 0.009)ftw 0. 
Discussion. There are a number of methodologies in 
current use to construct effective interactions for many- 
particle systems; among them, effective field theory 
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efficient. Non- unitary transformations might give faster 
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FIG. 3: (a) versus q for the lowest L — excited state 
of the A = 3 system. The inset shows Ji?' 9 - 1 versus q in a 
logarithmic scale. Symbols and units as in Fig. [2] (b) ' 
versus q for the L = ground state of the A = 4 system. 
The inset shows Ai?' 9 -' versus g for the new interaction in a 
logarithmic scale. The solid line is a linear fit to q = 1, 2, 3. 



(EFT) and the unitary-transformation method have a 
connection to the interactions discussed here. In EFT, 
the interaction is parameterized by contact terms (lead- 
ing order) and their derivatives. Our procedure to con- 
struct the g-renormalized contact interaction can thus 
be considered as leading-order EFT. Its poor conver- 
gence suggests that EFT treatments will require deriva- 
tive terms to accurately model trapped fermion systems. 

Our improved interaction has some connection with 
Suzuki's unitary regularization a method widely 

used in nuclear physics ID) 3 12 1- In Suzuki's approach, 
an effective interaction is determined by a unitary trans- 
formation of the Hamiltonian that decouples a subspace 
from its complementary subspace. In practice, the trans- 
formation is performed on the two-particle Hamiltonian, 
giving a transformed Hamiltonian that is block diagonal. 
This block diagonal structure guarantees that the energy 
eigenvalues are reproduced in the truncated subspace. 
Our effective interaction also reproduces the exact two- 
particle spectrum in a truncated subspace but has the 
advantage of being simple, i.e., separable. 

The unitary transformation of the two-particle Hamil- 
tonian cannot be carried out independently for all possi- 
ble pairings in the many-body Hamiltonian. When this 
transformation is applied to the many-particle system, it 
generates higher-order many-body interactions that are 
usually simply neglected. For our Hamiltonian, addi- 
tional correction terms would be required if we were to 
relate it to a unitary-transformed Hamiltonian. Rather 
than attempting to compute these correction terms, we 
have studied the convergence in the large g-limit, where 
our effective interaction coincides with the contact inter- 
action. By studying the convergence, one can assess the 
usefulness of many of the specific details of the different 
methodologies. For example, there are other choices of 
the many-particle space truncation that might be more 



convergence. The no-core-shell-model methodology 
is an example where a particular choice was made. 

Our method can be applied for interaction strengths 
away from unitarity, at the slight cost of inverting nu- 
merically a (q + l)-dimensional matrix. It may also be 
interesting to apply the method to uniform systems, us- 
ing the separability of the interaction in a plane-wave 
basis. One caveat is that we have only examined three- 
and four-particle systems. It will be important to confirm 
the exponential convergence when our interaction is used 
to calculate the spectra of systems with more particles. 
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